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ABSTRACT 

We present a sky simulatiorQ of the atomic HI emission line and the first ten ^^C^^O rotational 
emission lines of molecular gas in galaxies beyond the Milky Way. The simulated sky field has 
a comoving diameter of 500 Mpc, hence the actual field-of-view depends on the (user-defined) 
maximal redshift Zmax; e.g. for Zmax = 10, the field of view yields 4 x 4deg^. For all galaxies, 
we estimate the line fluxes, line profiles, and angular sizes of the HI and CO emission lines. The 
galaxy sample is complete for galaxies with cold hydrogen masses above 10® Mq. This sky simulation 
builds on a semi-analytic model of the cosmic evolution of galaxies in a A-cold dark matter (ACDM) 
cosmology. The evolving CDM-distribution was adopted from the Millennium Simulation, an A^- 
body CDM-simulation in a cubic box with a side length of 500 Mpc. This side length limits the 
coherence scale of our sky simulation: it is long enough to allow the extraction of the baryon acoustic 
oscillations (BAOs) in the galaxy power spectrum, yet the position and amplitude of the first acoustic 
peak will be imperfectly defined. This sky simulation is a tangible aid to the design and operation 
of future telescopes, such the SKA, the LMT, and ALMA. The results presented in this paper have 
been restricted to a graphical representation of the simulated sky and fundamental dN/ dz-analyzes 
for peak flux density limited and total flux limited surveys of HI and CO. A key prediction is that 
HI will be harder to detect at redshifts z > 2 than predicted by a no-evolution model. The future 
verification or falsification of this prediction will allow us to qualify the semi-analytic models. 
Subject headings: galaxies: high-redshift — galaxies: evolution — ISM: atoms — ISM: molecules — 
cosmology: theory 



1. INTRODUCTION 

The interstellar medium (ISM) is the bridge between 
the environment of galaxies and their newborn stars. 
Its atomic and molecular phases can be detected via 
emission lines. Typically studied lines include the HI- 
radio line (1.420 GHz rest-frame) and the rotational 
CO-lines in the (sub) millimeter spectrum (multiples of 
115.27 GHz). These lines characterize both the com- 
position and the dynamical state of the ISM, and their 
apparent frequency measures the redshift of the source. 
If the object sits at a cosmological distance, the red- 
shift is dominated by the expansion of the Universe and 
hence provides a distance measure. Therefore, observa- 
tions of HI and CO at high redshift are currently dis- 
cussed as a means of localizing high-redshif t galaxies 
(ICarilh fc Blaia[200l[Carini fc Rawlingi[200l . thus un- 
veiling an u nprecedented image of cosmic structure (e.g. 
lAbdalla et al. 2009). 

However, due to current sensitivity limitations, no HI- 
emission has yet been found beyond reds hift z = 0.25 
(jVerheiien et al.l2007t[Catinella et aO2008[ ). In contrast, 
CO-emission lines ha ve been detected in different sys- 
tems out to z « 6.4 () Walter et al.ll2004[ ). yet all these 
lines originate from atypical objects, such as ultra lu- 
minous infrared galaxies (ULIRGs) or quasi stellar ob- 
jects (QSOs), with the exception of two ordinary galax- 
ies recently detect ed in CO (2-l)-emission at z » 1.5 
(iDaddi et al.|[200l. Both H I ()Prochaska et al.|[2005l) and 
Hc> ([Noterdaeme et al.ll2008l ) have also been measured via 
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absorption of their respective Lyman lines against distant 
QSOs. Yet, the nature of the absorbing galaxies remains 
unknown. 

The discrepancy between the primordial astrophysical 
importance of cold gas in galaxies and its gravely lim- 
ited detectability at high redshift is a main driver for the 
design of many future radio and (sub)millimeter tele- 
scopes. Prominent examples are the Square Kilometre 
Array (SKA), the Large MiUimeter Telescope (LMT), 
and the Atacama Large Millimeter/submillimeter Array 
(ALMA), which are expected to detect HI and CO at 
high z. The optimization of these instruments and the 
planning of their surveys require robust predictions of 
the detectable signatures. Such predictions are available 
for t he continuum r adiation of a large sample of galaxies 
(e.g. lWilman et an [2b08). By contrast, high-redshift line 
emission of HI and CO h as only be en simulated for sin- 
gle galaxies (iBoomsma e t al. 200% iCombes et al.l 119991 : 
'Greve fc Sommer-Larsenll2008ll or simplistically extrap- 
olate d for a population of r egular galaxies (iBl ain ct al] 
&M ICarilh fc Blah] [200l lAbdalla fc R"OTdhrggl L2005). 
All these models ignore the implications of galaxy merg- 
ers, cooling flow suppression mechanisms, and other com- 
plex phenomena. Moreover, most line simulations exclu- 
sively consider either the atomic or the molecular gas 
phase of the ISM. This approach implicitly assumes that 
the other phase is negligible or that the masses in both 
phases evolve proportionally. Both of these assumptions 
contradict recent studi es of the co-evolution of HI an d 
H2 in regular galaxies (jObreschkow fc Rawlingsl[2009bf ). 
Finally, the line simulations cited above neglect cosmic 
large-scale structure. The time seems ripe for joint pre- 
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dictions of atomic and molecular emission lines in a sam- 
ple of galaxies large enough to probe cosmic structure. 

In this paper, we present a simulation of a sky field 
with a comoving diameter of 500/i^^Mpc. The actual 
field-of-view depends on the (user-defined) maximal red- 
shift Zmax; e.g. for Zmax = 1 the field of view yields 
~ 12 X 12deg^, or for Zmax — 10 the field of view 
yields ^ 4 x 4deg^. This simulation is obtained by con- 
structing a mock observing cone from a previously pre- 
sented galaxy simulation. The latter relies on the large- 
scale structure compute d by the Millennium Simula- 
tion (ppringel et al.ll2005D and an enhanced semi-analytic 
galaxy model (iCroton et "all 120061 : iDe Lucia &: BlaizotI 
l2007t iObreschkow et all 1200931) ! 

Section [2] explains the simulation methods. In Section 
[21 we provide a graphical illustration of the simulated sky 
field and extract dA^/dz-estimates for peak flux density 
limited surveys. Section |4] discusses some important lim- 
itations of the presented simulation. A non-exhaustive 
list of possible applications is provided in Section[5]along 
with a brief conclusion. The appendices show additional 
illustrations, list the parameters describing the analytic 
flts to the predicted dA^/dz-functions, and describe the 
on-line access to the simulation data. 

2. METHODS 

In this section, we describe the multiple simulation 
steps required to progress from a simulation of the evo- 
lution of cosmic structure to a static sky simulation. We 
have grouped this description into four steps, correspond- 
ing to four successive simulation steps. The first step 
(Section 12. 1|) contains all the simulation work presented 
in earlier studies. This work resulted in a catalog of 
^ 3 • 10^ evolving galaxies with detailed cold gas prop- 
erties. In the second step (Section 12.21) . this catalog is 
transformed into a mock observing cone, which repre- 
sents a virtual sky field. In the third step (Section 12. 3p . 
the intrinsic properties of the galaxies in this virtual sky 
field are converted into apparent line fluxes. In the fourth 
step (Sections 12.41 and 12. 5p , the line emission is reflned 
by the evaluation of line proflles and angular sizes of the 
line-emitting gas. 

2.1. Simulation of the ISM in '^3-10'^ evolving galaxies 

Here, we recapitulate the galaxy simulation presented 
in earlier studies. This simulation relies on three consec- 
utive layers: (i) a simulation of th e cosmic evolution of 
dark matter ([Springel et al.l [20051) : (ii) a semi-analytic 
simulation of t he evolution of galaxies on the dark 
matt er skeleton (jCroton et al.ll2006l : IDe Lucia fc BlaizotI 
I2007D : and (iii) a post-processing to split the cold hy- 
drogen masses associated with each galaxy into HI and 
H2 (Obreschkow et al. 2009a). 

For the dark matter sim ulation, w-e ado pted the Mil- 
lennium Simulation (|Springel et al.l 120051 ) . an A^-body 
dark matter simulation within the standard A-cold dark 
matter (ACDM) cosmology. This simulation uses a cu- 
bic simulation box with periodic boundary conditions 
and a comoving volume of {500 Mpc)^ . The Hub- 
ble constant was flxed to Hq = 100 ft.kms""'^ Mpc~^ with 
h = 0.73. The other cosmological parameters were cho- 
sen as fimattor = 0.25, ^Ibaryon = 0.045, f^A = 0.75, and 
as = 0.9. The simulation-box contains ~ 10^° particles 



with individual masses of 8.6 • 10* M©. This mass resolu- 
tion allows the identiflcation of str uctures as low in mas s 
as the Small Magellanic Cloud (see lSpringel et alll2005D . 

For the second simulation-layer, i.e. the cosmic evo- 
lution of the galaxies distributed on the dark mat- 
ter skeleton, we adopt the semi -analytic mo del of 
IDe Lucia fc BlaizotI (|2007[ ) (see also ICroton et al.i r2006). 
In this macroscopic model all galaxies are represented 
by a list of global properties, such as position, velocity, 
and total masses of gas, stars, and black holes. These 
properties are evolved using empirically or theoretically 
motivated formulae for mechanisms, such as gas cool- 
ing, reionization, star formation, gas heating by super- 
novae, starbursts, black hole accretion, black hole coales- 
cence, and the formation of stellar bulges via disk insta- 
bilities. The resulting virtual galaxy catalog (hereafter 
the "DeLucia-catalog" ) contains the positions, velocities, 
merger histories, and intrinsic properties of ~ 3 • 10^ 
galaxies at 64 cosmic time steps. The free parameters 
in the semi-analytic model were tu ned to various obser- 
vations in the local universe (see ICroton et all [2006). 
Therefore, despite the simplistic implementation and 
the possible incompleteness of this model, the DeLucia- 
catalog nonetheless provides a good flt to the joint lumi- 
nosity /colour /morphology_distribution of observed low- 
redshift galax i es JCole et al.l 120011 : iHuang et ahl 120031 : 
[Norberg et al.''2002^, the bulge-to-black hole mass re- 
lation (.Haring fc Rix 2004). the Tully-Fisher relation 
(jGiovanelli et al.lll997[ ). and the cold gas meta llicity as a 
function of stellar mass (jTremonti et al.|[200l . 

In this paper, we are particularly interested in the 
cold gas masses of the galaxies in the DeLucia-catalog. 
These cold gas masses are the net result of (i) gas accre- 
tion by cooling from a hot halo (dominant mode) and 
galaxy mergers, (ii) gas losses by star formation and 
feedback from supernovae, (iii) and cooling flow sup- 
pression by feedback from accreting black holes. The 
DeLucia-catalog does not distinguish between molecular 
and atomic cold gas, but simplistically treats all cold gas 
as a single phase. The atomic and molecular phases are 
therefore dealt with in the third simulation layer. 

The third simulation-layer, i.e. the subdivision of the 
cold hydrogen mass of each galaxy int o HI- and H2- 
distributions (IObreschkow et all l2009af l. relies on an 
analytic model for the mass-distributions of HI and 
H2 within regular galaxies. In this model, the column 
densities of HI and H2, Shi and Sjja respectively, are 
given by 



SHi(r) = 



SH.(r) = 



Sh exp(-r/rdisk) 
1 + ^moi exp(-1.6 r/rdisk) 
flmoi exp(-2.6r/rdisk) 
1 + ^moi exp(-1.6 r/rdisk) 



(1) 



(2) 



where r denotes the galactocentric radius, rdisk is a 
scale length, R'^^i is the H2/HI-mass ratio at the galaxy 
center, and Ejj is a normalization factor. We derived 
Eqs. ()l|2p based on a list of empirically supported as- 
sumptions, the most important of which are: (i) the 
cold gas of regul ar galaxies resides in a flat disk (see 
iLerov et"ari 120081 for local spiral galaxies. lYound 120021 
for local elliptical galaxies. iTacconi et al.ll2006l for galax- 
ies at higher redshifts); (ii) the surface density of the 
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total hydrogen component (HI+H2) is we ll described 
by an axially symmetric exponential profile ()Lerov et alj 
[2001 : (iii) the local H2/HI -mass ratio scales as a power 
of the gas pressure of the ISM outside molecular clouds 
(jBlitz fc Rosolowskvl[2006l ). 

Using Eqs. (|l|2p . we can characterize the HI- and 
H2-content of every simulated galaxy in the DeLucia- 
catalog. The resulting hydrogen simulation successfully 
reproduces many local observations of HI and H2, such 
as mass functions (MFs), mass diameter relations, and 
mass-velocity relations (Obrcschkow ct al. 2009a). This 
success is quite surprising, since our model for HI and 
H2 only introduced one additional free parameter to 
match the observed a .verage space density of cold gas 
in the local Universe (jObreschkow et aLll2009al ). A key 
prediction of this simulation is that the H2/HI-ratio of 
most regular galaxies increases dramatically with red- 
shift, hence causin g a clear signature of cosmic dow nsiz- 
ing in the H2-MF (jObreschkow fc Rawlindl2009bD . 

Despite its consistency with existing observations, 
we emphasize that the presented model for the cos- 
mic evolution of HI and H2 is simplistic and uncer- 
tain. In particular at high 0, the model assumptions 
may significantly differ from the reality. For exam- 
ple, high-z galaxies are likely to be more disturbed due 
to higher merger rates and long dynamic time scales 
compared to their age. There is also evidence that 
cold gas disks become more turbulent with redshi ft (e.g. 
IForster Schreiber et al.l [20061: IGenzel et all [2001 . Un- 
certainties from these and other model limitations are 
di scussed briefly in Section jH and in depth in Section 6 
of jObreschkow et al.l (|2009ar 

2.2. Building a mock observing cone 

We shall now describe how the cubic simulation-box is 
transformed into a virtual sky field. This procedure can 
be regarded as a fourth simulation-layer on top of the 
hierarchical simulation described in Section [2. II 

The met hod adopted h ere cl osely follows the one de- 
scribed bv iBlaizot et al.l (|2005D . namely the building of 
a chain of replicated simulation-boxes along the line-of- 
sight, as shown in Fig. [TJ At any position in this chain, 
the galaxies are drawn from the cosmic time in the simu- 
lation, which corresponds to the look-back time, at which 
the galaxy is seen by the observer O. Since our galaxy 
simulation uses 64 discrete time steps, we describe each 
galaxy in the cone by its properties at the closest avail- 
able time ste]:[3, in terms of redshift. This defines the 
spherical shells of identical cosmic time, which are sep- 
arated by dashed lines in Fig. [T] The relatively narrow 
redshift-separation of these shells ensures that the as- 
signed galaxy properties cannot differ significantly from 
the properties at their exact look-back time. 

The same galaxy can appear once in every box in 
Fig. [1] but with different intrinsic properties due to 
the cosmic evolution. However, the position of the re- 
peated galaxy in comoving coordinates will be very sim- 
ilar, which can result in s purious radial featu res for the 
observer O (see Fig. 1 in IBlaizot et al.l [20051) . To sup- 
press this effect, we randomize the galaxy positions by 

^ The galaxy properties cannot readily be interpolated between 
two successive time steps, since a galaxy at any time step may have 
several progenitors. 
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Fig. 1. — Schematic illustration of the construction of a mock ob- 
serving cone (shaded region) from a chain of replicated simulation 
boxes (solid square boxes). The galaxies in the replicated boxes 
are drawn from the discrete cosmic-time step of the simulation, 
which best corresponds to their cosmological redshift relative to 
the observer O. Galaxies from the same discrete time step there- 
fore lie within spherical shells around the observer (indicated by 
the dashed-lines). 



applying random symmetry operations to each box in 
the chain. These operations consist of 90 deg-rotations, 
inversions, and continuous translation^ Applying these 
symmetry operations also removes the non-physical peri- 
odicity of 500 hr^ Mpc associated with the side length of 
the periodic simulation-box. But we emphasize that ap- 
plying the symmetry operations does not provide infor- 
mation on scales larger than the simulation-box. Symme- 
try operations can, however, introduce unwanted small- 
scale density variations at the interface of two neighbor- 
ing boxes. Thes e and other hmitati ons of this method 
are discussed by IBlaizot et al.l ()2005[) . 

From the randomized chain of replicated simulation- 
boxes, an observing cone can be extracted (shaded re- 
gion in Fig. [TJ. Each galaxy in this cone is projected 
onto the celestial sphere centered about the vernal point 
(RA = 0, Dec = 0). The Euclidian projection formulas 
for arbitrary large angles are 



RA — arctan 



DEC = arctan 



(3) 



(4) 



where Tx , ry , and rz are the comoving coordinates of the 
galaxy relative to the observer (see Fig. [IJ. The "cos- 
mological redshift" z of each galaxy is computed directly 
from its comoving distance Dc — (^x + '"y + ^z)^^^: while 
the Doppler-shift corrected "apparent redshift" is com- 
puted as z = z -I- V,./c, where is the peculiar recession 
velocity of the galaxy relative to the Hubble flow. 

Fig.[l]shows that the opening angle of the virtual sky 
field is set by the maximal comoving distance I?c,max via 



2 arcsin ■ 



(5) 



where Sbox is the comoving side length of the simulation- 
box. Given a value of Sbox and a choice of cosmological 
parameters, Eq. ([5]) implies a one-to-one relation between 
if and the maximal redshift 

Fig. [2] displays the relation between ip, Dc.max, and 
•^max for the cosmological parameters of the Millennium 
Simulation (Section 12. ip and three different side lengths 

^ Translations can be applied because of the periodic bound- 
ary conditions imposed on the simulation box of the Millennium 
Simulation. 
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Sbox- The choice Sbox = 500/i~^Mpc (sohd line) corre- 
sponds to the box of the Millennium Simulation. Sbox = 
62.5 Mpc (dashed line) corresponds to the "Milli- 
Millennium" Simulation, a small test version of the Mil- 
lennium Simulation. Sbox — ^h^^ Gpc (dash-dotted line) 
corresponds to the giant simulation-box of the Horizon- 
47r Simulation, a dark matter stimulation with ~ 10- 
times l ess mass resolution than the Millenniu m Simu- 
lation UPrunet et al.ll2008t iTbvssier et all 120081: see also 
Section US. 
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Fig. 2. — Relation between the maximal comoving distance 
^C,max or maximal redshift Zmax and the opening angle tp of 
the mock observing cone. The different lines correspond to 
the box sizes of the Millennium Simulation (solid), the Milli- 
Millennium Simulation (dashed), and the Horizon-47r Simulation 
(dash-dotted). 



The Hl-line or "21 cm- line" at a rest- frame frequency 
of 1.420 GHz originates from the photon- mediated tran- 
sition between the two spin states of the proton- 
electron system in the electronic ground state. The 
upper spin state has a low spontaneous decay rate of 
/ = 2.9 • 10~^^ Hz. This frequency is about 5 or- 
ders of magnitude smaller than that of HI-HI collisions 
(jBinnev fc Merrifieldlll998() . Hence, the two spin states 
are in thermal equilibrium with the kinetic state of the 
gas, which implies a spin-temperature far above the spin 
excitation temperature T^x ~ 0.07 if. Therefore, the 
spin systems are in the high-temperature limit, where 
3/4 of all systems are in the upper (three-fold degener- 
ate) state. The radiative power emitted per atom can 
therefore be calculated as Lhi = 0.75 f hp 1.4 GHz {hp is 
the Planck constant), which readily reduces to Eq. 

We have neglected Hl-self absorption, since this seems 
to affect only r nassive spiral g alaxies wh en observed 
almost edge-on ()Ferrierd l200lt IWalll [20061 ). This as- 
sumption should also be valid for high- redshift galaxies, 
since their Hl-masses were not much larger than today, 
as can be inferred from Lyman-a absorption measure- 
ments against distant quasars (Lah ct al. 2 00 7|) and as 
is pred icted by our simulation (jObreschkow fc Rawlingsl 
l2009bD . 

HI in collapsed structures, i.e. in galaxies, is generally 
warm enough (> 50 K) that the Cosmic Microwave Back- 
ground (CMB) can be safely neglected as an observing 
background for all galaxies at z < 10. Only in the in- 
tergalactic medium (IGM) during the cosmic Epoch of 
Reionization (EoR) c an HI appear in 21 cm-absorption 
against the CMB (e.g. llliev et aLll2002[ ). 



2.3. Assigning apparent line fluxes 
We shall now assign apparent line fluxes to each galaxy 



in the mock observing cone constructed in Section [2?2] In 
general, the frequency-integrated line flux S of any emis- 
sion line, can be computed from the frequency-integrated 
luminosity (= intrinsic power) via 



S ■ 



(6) 



where Z?l — (1 + z) Dc is the luminosity distance to 
the source. We note that Eq. ^ takes a different form 
for velocity-integrated fluxes and brightness tempera- 
tures. For reference, a summary of all relations between 
frequency-integrated and velocity-integrated fluxes, lu- 
minosities, and brightness temp eratures has b een com- 
piled in Appendix A of Obresch kow et all (|200 9b). 

In Scctions l2.3.1l and l2.3.21 we shall now summarize the 
models used to estimate the frequency-integrated lumi- 
nosities L (hereafter simply "luminosities") of different 
emission lines. 

2.3.1. Conversion of Hl-mass into Hl-line luminosities 

We evaluate the Hl-luminosities Lhi from the HI- 
masses Mhi of the simulated gala xies via the standard 
conversion (e.g. iMever et al]|2004D 



6.27-10" 



HI 



M^ 



(7) 



2.3.2. Conversion of H2 -mass into CO-line luminosities 

We derive the CO-line luminosities Leo from the H2- 
masses of the simulated galaxies. We only consider the 
radiation emitted by the most abundant CO-isotopomer, 
12(^16q^ when relaxing from one of the rotational states 
J = 1, 10 to the state J — 1. The radiation frequency 
of such a decay is J ■ vco, where vco = 115.27 GHz is 
the rest-frame frequency of the transition J = 1— )-0. 

The conversion between H2-masses and CO- 
luminosities is a highly nuanced affair with a long 
history in millimeter astronomy. We therefore pre- 
sented an in-dep t h ana lysis of this conversion in 
lObreschkow et al] (|2009b( ) and introduced a model to 
estimate the different luminosities Leo, J of the galaxies 
in our simulation. This model respects the following 
physical mechanisms: (i) molecular gas is heated by 
the CMB, starbursts (SBs), and active galactic nuclei 
( AGNs) ; (ii) molecular clouds in dense or inclined galax- 
ies can overlap; (iii) very dense gas is smooth instead of 
clumpy; (iv) the metallicity varies amongst galaxies and 
changes with redshift; (v) CO-luminosities are always 
detected against the CMB. We shall apply this model 
in the present paper. Limitations and uncertainties are 
discussed in Section |4?3l 



2.4. Emission line profiles 

Having assigned an integrated line flux to each galaxy 
in the simulation, we can now refine their attributes, 
by characterizing each line with a profile. To this end 
we depart from the edge-on line profiles evaluated in 
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lObreschkow et all (|2009aD for each galaxy in the simu- 
lation. We represented those profiles by normalized flux 
densities VE'(T^), where V is the velocity measured in the 
rest-frame of the center of the observed galaxy. The nor- 
malization condition, 'i'{V)dV = 1, implies that '^/{V) 
only needs to be multiplied by the velocity-integrated 
flux (in units of Jykm/s) in order to obtain actual flux 
densities (in units of Jy) . For each galaxy we calculated 
two profiles 5'(^), one for the Hl-component and one for 
the H2-component (associated with C0j3- This calcula- 
tion relied on a detailed mass model based on the halo, 
disk, and bulge of the galaxies, combined with our model 
for the HI- and H2 -surface densities given in Eqs. (|l|2p . 
For practical reasons, the resulting line profiles were re- 
duced to five parameters (see Fig. |S]) : the normalized flux 
density at the line center ^'o; the normalized peak flux 
density ^max, usually corresponding to the two peaks of a 
double-horn profile; the line width Wpcak between the two 
peaks of the double-horn profile; the line width W50 at 
the 50-percentile level of the peak fiux density; and the 
line width W2q at the 20-percentile level. The original 
normalized line profile can be approximately recovered 
from these parameters using the formulas in Appendix 

El . . 

The remaining task consists of correcting the line pro- 
files for the inclination i of each galaxjQ. i is defined as 
the angle between the line-of-sight and the galaxy's rota- 
tion axis; hence i = deg corresponds to face-on galaxies 
and i = 90 deg corresponds to edge-on galaxies. In the 
absence of a random gas velocity dispersion, apparent 
line widths w°^^ could be computed from the edge-on 
line widths w, via w°^^ = w ■ smi, while apparent nor- 
malized flux densities would scale as = ^P/ sini. 

Here, we assume that the cold gas has a random, 
isotropic velocity dispersion characterized by a Gaus- 
sian velocity distribution in each space dimension. The 
observed line profile of a face-on galaxy (j = deg) 
then takes the shape of a Gaussian function. Under 
no inclination can the line profile become more nar- 
row than this Gaussian function, nor can the normal- 
ized line fiux densities become higher than the peak of 
this Gaussian function. Let a be the standard-deviation 
of the Gaussian velocity dispersion. Then the minimum 
line widths are given by u;^^'' = 2y/-2 ln(0.2) « 3.6 tr, 
= 2^-2 ln(0.5) « 2.4 cr, and = 0> and the 

maximum normalized fiux densities are — {<j^/2t:)~^ 
and ^-max = (cr\/27r)~^ 

In addition, the maximal normalized fiux density ^max 
cannot differ from the central normalized flux density '^q 
by an arbitrarily large amount, due to the line profile 
smoothing. Explicitly, the slope in the emission line be- 
tween the points ^Pq and ^max cannot exceed the max- 
imal slope of the Gaussian velocity function, which is 
equal to 0.24 cr~^. This requirement translates into an 
upper bound for ^-^ax equal to ^'^'^ + 0.12 ct"^ w°oik- ^ 
set of equations respecting all of these conditions is given 



^ In lObreschkow et al.l l l2009al ). we called these two profiles 'J'jji 
and "l/cGi respectively. 

* The simulated DeLucia-catalog does not provide galaxy orien- 
tations. We therefore assign inclinations randomly between Odeg 
(face-on) and 90 deg (edge-on) according to a sine-distribution. 



by 

W2o' = (w20-3.6cr)-sini + 3.6cr, (8) 
t(;°b'' = (w50 - 2.4 cr) • sin i-F 2.4 cr, (9) 

^«poak = ^«pcak-sini, (10) 

<^=minfi^,^), (11) 
\sim crv 27r/ 

^S.^x-minf^,^,*--fH!!^y (12) 
\ smi cr V27r cr I 

For all the line profiles, we here adop t cr = 8kms~^ 
to remain consistent with Ob reschkow et al.l (|200 9aV 
We note, however, that high -redshift galaxies may 
have higher velocity dispersions (jForster Schreiber et al.l 
I2006f ) perhaps due to an intense ongoing accretion of gas. 
Another limitation of the presented line model is that all 
CO-lines have by definition the same line shape. This 
assumption nevertheless approximately agrees with si- 
multaneous observations of differen t emission lines (e.g. 
iWeiss et al.l[2007l: IGreve et al.l[2009l) . 

2.5. Angular sizes 

We shall finalize our simulation of line-emitting galax- 
ies by ascribing an angular distribution (flux per unit 
solid angle) in the sky to each emission line. To this 
end we assume that, for a fixed galaxy, the angular dis- 
tributions of the Hl-flux and the different CO-fluxes are 
proportional to the angular distributions of HI and H2, 
respectively. We emphasize, however, that the conver- 
sion factor between H2-masses and CO-fiuxes depend on 
the J-level of the CO-transition and on the galaxy, such 
as outlined in Section r2. 3. 21 

Given this assumption, we can infer the line fiux densi- 
ties per unit solid angle from the surface densities I]Hi(f ) 
and Tj-n^ir) given in Eqs. (|ll2p . We only need to normal- 
ize these densities to the respective line fluxes and to 
replace the scale radius rdisk by the apparent scale ra- 
dius rdisk/-DA, where Z?a = (1 + z)~'^ Dq is the angular 
diameter distance. 

The shapes of the line-emitting regions need to be cor- 
rected for the inclinations i. If <7hi,o and qn^fi respec- 
tively denote the intrinsic axis ratios of the atomic and 
molecular gas in galaxi es, then the ap parent axis ratios 
are given by (Eq. (1) in lKannappan e t al. 2002.) 

2 2 • I 2 • 2 • /I o^ 

gni^cos I + (7hi,o sm (13) 

gg^=cos2i + g2^_^ sin^i. (14) 

These relations satisfy qm = = 1 , if « = deg (face- 
on), and qm = qm,o and q^^ = gHs.o, if « = 90 deg 
(edge-on). 

We assume that all galaxies have the same values for re- 
spectively qm,o and gH2,o- In local spiral galaxies, we typ- 
ically find qm.o — 0.1, as can be seen from high-resolution 
maps of the edgc-on spiral galaxies NGC 891 and NGC 
4565 (Rupcn 1991). To our knowledge, no reliable es- 
timate of (7h2,o for disk galaxies is available. However, 
simultaneous 00(1-0) and optical observations revealed 
that the density of stars in nearby galaxies correlates 
with the density of molecular gas, when averaged over 
sufficiently large (> kcp) scales fiRichmond fc Knapnl 
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Fig. 3. — Longitudinal cut through the simulated mock observing cone. The cut slice has a thickness of lOMpc and is represented in 
comoving coordinates. For illustration purposes, the slice has been wrapped in four parts, which variably read from the left to the right 
and vice versa. The dots represent gas-rich galaxies and the coloring shows their H2/HI-ratio, from (blue) to 2 (red). 



[I986t iLerov et alj [20081 ) . probably as a natural conse- 
quence of the formation of stars from molecular gas. 
Therefore, we shall assume that the intrinsic aspect ra- 
tio of molecular gas qH2,o is identical to that of stel- 
lar disks. The latter is ~ 0.1, as can be seen in the 
sample of 34 nearb y edge-on spiral galaxies studied by 
iKregel et all (|2002D . We therefore adopt qhs.o =0.1. 
We stress that qm,o — 'i'Ha.o does not contradict the 



fact that characteristic scale radii and scale heights of 
Hl-distributio ns are generally l arger than those of H2- 
distributions ()Lerov et al.ll2008f ). 

The assumption of constant values for qm,Q and gH2.o 
for all galaxies at all redshifts may not be verified at 
high redshifts due to time-scale arguments. In fact, a 
significant fraction of the galaxies at 2; > 5 may have 
an age comparable to their dynamical time-scale. Their 
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Fig. 4. — Illustration of the galaxies in the redshift range z = 1.0 — 1.1 in a small field of 1 arcmin-^. The full field of view of the observing 
cone is 60,000-times larger than this example. The gradual coloring represents integrated line fluxes per unit solid angle for HI (left) and 



CO(l— 0) (right). The different color tones for CO represent the brightness temperature intensity ratio ^cO(5-4)/^CO 



)(i-0)- 



The white 



contours around Hl-sources represent iso-density curves of CO at the 50-percentile level of the full CO density scale and vise versa. 



cold gas distribution might therefore be thicker and less 
circularly symmetric than the flat gas disks seen today. 
However, no reliable estimates of qHi.o and <7h2,o beyond 
the local Universe are available today. 



3. RESULTS 

The simulated mock observing cone can be accessed 
on-line as described in Appendix [Bl This section starts 
with a graphical illustration of the simulated mock ob- 
serving cone. As a second step, we analyze the pre- 
dicted number of galaxies detected in an idealized line- 
survey with a constant peak flux density limit. Spe- 
cific predictions for particular surveys with radio and 
(sub)millimeter telescopes, such as the SKA, the LMT, 
and ALMA, shall follow in forthcoming studies. 

3.1. Graphical overview 

By successively applying the simulation steps described 
in Section [5J we have constructed an observing cone of 
line emitting galaxies. Fig. [3] shows a longitudinal slice 
of this cone with a thickness of 10 Mpc. This slice corre- 
sponds to a diagonal cut, as illustrated in Fig. [3] and has 
an opening angle of 5.8 deg. Each pixel inside this slice 
corresponds to a galaxy. The structure of the cosmic web 
appears clearly, as well as the increasing filamentarity of 
this structure with decreasing redshift. The color scales 
represents the H2/HI-mass ratios of the galaxies. We 
can clearly recognize the pressure-driven cosmic decline 
of this ratio (see Obreschkow & Rawlings 2009^- 

The mock observing cone translates into a virtual sky 
field when projected onto a sphere using Eqs. pi4p . Fig.|l] 
displays the HI- and CO-flux densities of galaxies be- 
tween z = 1 and z = 1.1 in a small extract of this 
virtual sky. The surface densities of the galaxies have 



8 



Obreschkow et ah 



been modeled using Eqs. (jll2p in the manner described 
in Section 12.51 The more massive galaxies in the field 
reveal ring-like Hl-distributions with CO-rich central re- 
gions. By contrast, some of the smaller galaxies with 
low surface brightness (LSB), have most of their HI in 
the center, with nearly no detectable CO. In general, 
CO (and hence H2) is more compact than HI. All these 
simulated features compare well to observed radial distri- 
butions (averaged on suitably large s cales) of HI and CO 
in nearby galaxies fLer oy et al.l l2008') as demonstrated in 
[obresch kow et al. (2009a|) andlObreschkow & Rawlings 
(|2009af) . Nevertheless, we emphasize that the idealized 
nature (e.g. axial symmetry) of the simulated galaxies 
misses out a list of structural features found in many real 
galaxies, such as spiral arms, warps, and lopsidedness. 

The coloring of the CO-surface densities in Fig. |4] 
(right) represents the CO(5~4)/CO(1-0) line ratio in 
terms of brightness temperatures. For normal galax- 
ies, without a particular source of heating, this ratio is 
much smaller than unity (yellow coloring), however for 
some galaxies with strong heating by an ongoing SB or 
AGN, the higher order lines can get excited (red color- 
ing). These me chanisms and our model t o assess them 
are discussed in [Obreschkow et all (|2009b[ ). 

Fig. iniin Appendix [Cl shows a 3-times larger sky field 
than Fig. U at the three redshifts 2; w 1, z w 3, and z w 6. 
The progression from z « 1 to 2; « 6 in Fig. |9] reveals 
two notable features. Firstly, galaxy sizes decrease with 
redshift. In fact, the angular diameter distances at z 1 
and z w 3 are virtually identical according to the cosmol- 
ogy adopted in this paper (Section [2.ip . Therefore, the 
galaxy sizes of these two virtual sky maps can be com- 
pared directly. The angular diameter distance at z w 6 
is about 25% smaller, hence the same physical scales ap- 
pear slightly oversized. The size evolution of the galaxies 
reflects the cosmic evolut ion of the volume/m ass-ratio of 
the dark matter haloes (jCunn fc Gottill972D . We dis- 
cussed the impact o f this evolution on the surface de nsi- 
ties of HI and H2 in I Obreschkow fc Rawlingsl (|2009aD . 

Secondly, the CO(5-4)/CO(l-0) line ratios of CO- 
detectable galaxies increases with redshift. In fact, at 
z w 6 no galaxy with a line ratio significantly below unity 
(i.e. with yellow coloring) can be seen. This feature re- 
lies partially on the compactness of t he galaxies, which, 
according to our model for CO-lines ([Obreschkow et all 
l200 9bf) , allows an efficient heating by star formation. An 
additional reason for the absence of low CO(5-4)/CO(l- 
0) line ratios at z > 6 is that molecular gas in galaxies 
with no significant star formation and no AGN will be 
hardly detectable in GO due to its near thermal equilib- 
rium with the CMB. 

3.2. dN/dz for a peak flux density limited survey 

In the simulated observing cone, we can readily count 
the number of galaxies per redshift interval with line 
fluxes above a certain threshold. This diV/dz-analysis 
is a key step towards a prediction of the number of 
sources detectable with any particular telescope and sur- 
vey strategy. In this section, we focus on the number of 
sources detected in a peak flux density limited survey, 
and we restrict the presented results to the HI, CO(l-O), 
CO(5-4), and CO(10-9) emission lines. Results for other 
CO emission lines and/or for integrated fiux limited sur- 
veys are presented in Appendix ID] 



Fig.[S]shows the dA^/dz-functions for six different peak 
flux density limits, logarithmically spaced between 1 mJy 
and 10 njy. Peak flux densities for each source and emis- 
sion line are calculated in the way described in Section 
12.41 This method accounts for the different gas distri- 
butions, rotation curves, and inclinations of the galax- 
ies. Every source with a peak fiux density above the 
peak flux density limit is considered as detected, while 
all other sources are considered as non-detected. Differ- 
ent aspects of Fig. [5] will be discussed in detail over the 
following paragraphs. 

3.2.1. Cosmic variance 

The simulated dA^/dz-functions shown in Fig. [S] (solid 
and dotted lines) correspond to a bin size of Az = 0.1 
and a sky field of 4 x 4 deg^, extracted from one partic- 
ular realization of the mock observing cone, that is one 
random choice of symmetry operations for the replicated 
simulation boxes (see Section [2^ . The wiggles visible 
in the simulated dA/dz-functions are physical. Similar 
wiggles can be expected for a real sky survey of a sky 
field of 4 X 4 deg^ with a redshift bin size of Az = 0.1. 
The fact that the amplitude of those wiggles does not 
decrease as 1/ dN/dz clearly uncovers the presence of 
the cosmic large-scale structure, also visible in Fig. |3l 

To quantify the effects of cosmic variance, we now con- 
sider the dA^/dz-functions extracted from five different 
random realizations of the mock observing cone. Fig. [6] 
shows the corresponding dA^/dz-functions for a peak flux 
limited Hl-survey with a peak fiux limit of 1 /iJy. Each 
function uses a bin size of Az = 0.1 and a small sky 
field of 1 X 1 deg^ in order to make the effects of cosmic 
structure obvious. As a rough estimate the log-scatter 
between the different dA'/dz-functions is about 0.1 dex. 
From this small scatter we conclude that cosmic variance 
is, in most cases, negligible compared to the uncertainties 
of the semi-analytic galaxy model. 

However, the comoving volume per unit solid angle and 
unit redshift varies as a function of redshift. Therefore, 
the scatter due to cosmic variance varies with redshift. 
It is largest at the lowest redshifts (z < 0.5), where the 
comoving surface per unit solid angle is small, and at 
the highest redshifts (z > 5), where the radial comov- 
ing distance per unit of redshift is small. In these red- 
shift regimes effects of cosmic variance should therefore 
be estimated, when comparing simulated data to obser- 
vations. 

3.2.2. Completeness 

The semi- analytic simulation is complete for galaxies 
with total hydrogen masses (HI+H2) above ^ 10^ M©. 
Galaxies with smaller hydrogen masses typically sit at 
the centers of halos with less than 20 particles in the 
Millennium Simulation, which cannot be reliably identi- 
fied. Therefore, the number of these low-mass galaxies is 
heavily underestimated by the simulation. In principle, 
the sensitivity of a line survey can be high enough that 
galaxies in this incomplete region of the simulated MF 
are detected. 

In fact, for each emission line and each sensitivity limit, 
there is a critical redshift Zc, below which sources of 
the incomplete part of the simulated galaxy MF will be 
detected at a sufficiently high rate that the number of 
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Fig. 5. — dA^/dz-plots for the emission lines of HI, CO(l-O), CO(5-4), and CO(10-9) at different peak flux density limits. Solid lines 
represent the simulated data in the regime where the simulated galaxy sample is complete. Dotted lines represent the cases where the flux 
density limit is low enough to detect sources in the incomplete parts of the galaxy mass function (see Section 13.21 1. These dotted lines 
should be considered as lower limits. To increase the readability, the transitions from the complete regime (solid lines) to the incomplete 
regime (dashed lines) have been marked with stars. All simulated functions (solid and dotted lines) use a redshift bin size of Az = 0.1 and 
a sky field of 4 X 4 deg^, corresponding to a maximal redshift 2max = 10. Dashed lines represent analytic fits of Eq. I I15I I. The respective 
fitting parameters are listed in table[T] Colors correspond to the peak flux density limits shown in the panel for CO(l-O). 
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Fig. 6. — Effects of cosmic variance on a peak flux limited Hl- 
survey with a flux limit of 1 /.tJy. The five lines show the dN/dz- 
functions extracted from five dist inct random realizations of the 
mock observing cone (see Section I2.2|l . Each function uses a bin 
size of Az = 0.1 and a small sky field of 1 X 1 deg^ in order to 
illustrate the effects of cosmic variance. 

real detections will significantly exceed the number of 
simulated detections. Tlierefore, the simulated dN/dz- 
functions at z < Zc must be considered as lower limits. 



We shall here define Zc as the redshift, below which the 
fraction of sources with hydrogen masses below 10* Mq 
becomes larger than 1%. This may seem to be a very 
low threshold, but we stress that the incomplete part of 
the simulated hydrogen MF misses a significant fraction 
of the "real" sources and hence the value of dN/dz may 
be underestimated by much more than 1% if z < z^. 
In Fig. [5] the values of Zc have been marked with stars, 
and the simulated d A^/dz-functions at z < Zc have been 
plotted as dotted lines instead of solid ones. 

3.2.3. Parametrization of the dN/dz-plots 

The simulated diV/dz-functions can easily be recov- 
ered from the on-line database of the sky simulation (see 
Appendix [B]) . Alternatively, we also approximated the 
dN/ dz-functions by analytic fits of the form 

= 10=^ • z^= • exp(-C3 • z), (15) 

deg"^ 

where ci, C2, and C3 are free parameters. The best param- 
eters in terms of an rms-minimization are shown in Table 
[1] for various emission lines detected with different limits 
for the peak flux densities and integrated fluxes. An- 
alytic dA^/dz-functions for intermediate flux limits can 
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be approximately inferred by linearly interpolating the 
parameters ci, C2, and C3. 

3.2.4. Basic conclusions 

An important conclusion from Fig. [5] is that Hl-surveys 
at high redshift (z > 2) wih be difficuh compared to 
CO-surveys. For example, in order to detect the same 
number of sources at 2: sa 4, an Hl-survey will need to be 
approximately 10-times more sensitive than a CO(l-O)- 
survey and over 100-times more sensitive than a C0(5- 
4)-survey. However, we emphasize that with single dish 
receivers this effect is partiahy offset by the fact the at 
the instantaneous field of view increases as , where A is 
the observed wavelength. This scaling implies that single 
dish Hl-surveys are likely to yield much more integra- 
tion time per unit solid angle, which effectively increases 
the sensitivity. Furthermore, new techno logies, such as 
aperture arrays (jCarilh fc R,awlingsl [20?)l. are currently 
being developed to provide an enormous instantaneous 
field of view and hence a very high effective sensitivity 
for the detection of HI. 

Fig. [7] (left) shows a comparison of the simulated 
dA^/dz-functions for different emission lines observed 
with an identical peak flux density limit of 1 /iJy. The 
flat slope of the dA^/dz-function for CO(10-9) reflects 
that this line is boosted by SBs, which were more abun- 
dant and effective (more c ompact galaxies) at high z (see 
lObreschkow et al.ll2009H ). CO(l-O) reveals the steep- 
est slope of all the CO-lines in the dA^/dz-plot. On one 
hand, this feature indicates that local galaxies are domi- 
nated by low-order excitatio ns of the CO-mole cule, con- 
sistent with empirical data (jBraine eFallflQQl . On the 
other hand, CO(l-O) becomes nearly invisible in normal 
(i.e. no AGN, no SB) galaxies at high redshift (z > 7) 
due to a near thermal equilibrium between the molec- 
ular gas_and_the CMB (lO breschkow et al.| '"2009b; see 
also lCombes et al.lll999llPap adopoulos et al. 2000). The 
even steeper slope of the dA^/dz-function for HI origi- 
nates from the c osmic decline of the H2 /H l -ratio i n galax- 
ies described in lObreschkow fc Rawlingsl ()2009bf ). 

Fig. [7] (right) shows the same dA^/dz-functions as Fig. [7] 
(left), but for the case of no galaxy evolution. These 
functions were obtained by constructing a mock observ- 
ing cone using only the simulation-box at z = 0. The 
comparison of Fig. [7] (left) to Fig. [7] (right) reveals that 
HI at high redshifts will be much harder to detect than 
predicted by a no-evolution model. Qualitatively, the 
same conclusion applies to low-order CO-emission lines, 
but the effect is less significant. In contrast, our sim- 
ulation predicts that the high-order CO-emission lines 
will be easier to detect than suggested by a no-evolution 
model, since these lines will be strongly boosted by SBs 
at high redshift. 

4. DISCUSSION 

4.1. Limitations of the galaxy simulation 

Our simulation is inevitably bound to the ACDM cos- 
mology with the cosmological parameters given in Sec- 
tion [21 The empirical uncertainty of these parameters 
may be a source of systematic errors in our predictions. 
To analyze the errors associated with the uncertainty of 
the Hubble constant, we can study the change of our 



predictions in the linear expansior0 of h. This analysis 
shows that varying h between 0.6 and 0.8 does not sig- 
nificantly affect the dA^/ dz-functions , i.e. n ot more than 
a factor 2. Additionallv. iWang et al.l (|2008[ ) showed that 
the lower value for the fiuctuation amplitude erg found 
by WMAP-3 compared to the value used in the Millen- 
nium simulation is almost entirely compensated by an 
increase in halo bias. Caution should nevertheless be 
applied when relying on predictions from a single cosmo- 
logical model. 

An additional limitation of the Millennium Simula- 
tion is the mass resolution of 8.6 • 10^ M0 per particle. 
This mass scale sets the completeness limit in our hy- 
drogen simulation to -M"hi + -^Ha ~ 10* M0 (Section 
13.21 lObreschkow et al.ll2009aD . Moreover, galaxies with 
Mhi -I- ^ 10^ M0 normally sit a the centers of 

dark matter halos with poorly resolved merger histories. 
Therefore, their properties m ay not have converged in 
the semi-analytic simulation ([Obreschkow et aD l2009a[ : 
ICroton et al1[2006l ). 

A long list of limitations associated with the semi- 
analytic galaxy simulation and our post-processing to 
assign e xtended HI- and H 2 -prope rties has been consid- 
ered bv lObreschkow et al.l ([2009a[ ). The bottom line of 
this discussion is that, at z > 5, the simulation becomes 
very uncertain because the geometries and matter con- 
tent of regular galaxies are virtually unconstrained from 
an empirical viewpoint. The young age and short merger 
intervals of these galaxies may, in fact, have caused 
them to deviate substantially from the simplistic disk- 
gas model. At z < 5, the predictions of our HI- and H2- 
properties are more certain, as they are consistent with 
available observations. For example, two measurements 
of C0(2-l)-line em ission in regular galaxies at z « 1.5 
(iDaddi et al.ll2008D are consistent with th e H2-MF at this 
redshift ([Obreschkow fc Rawlingsl [2009bD . Furthermore, 
the predicted comoving space density of H2 evolves pro- 
portional ly to the obser ved space density of star forma- 
tion (e.g. lHopkinsll2007[ ) within a factor 2 out to at least 
z = 3. At z = 0, the simulated HI-MF and CO(l-O)- 
lu minosity function a re co nsistent with the o bservations 
of lZwaan etall ([2005D and lKeres et~all ([2003D . Addition- 
ally, the local sizes an d line widths of HI and C O match 
the local observations ([Obreschkow et al1l2009aL and ref- 
erences therein). 

We shall now highlight some specific limitations asso- 
ciated with the emission lines considered in this paper. 

4.2. Limitations specific to the Hl-line 

We emphasize that at high redshift, the simulated cos- 
mic Hl-space density 57hi falls below the inferences from 
Lyman-a absorption against distant QSOs by a factor 
^ 2. As mentioned in lObreschkow fc Rawlingg ((2009b) . 
this could refiect a serious limitations of the semi-analytic 
models implied by the treatment of all cold hydrogen 
(HI4-H2) as a single phase. Consequently, our dA/dz- 
predictions for HI could be slightly pessimistic. If we 
believe the empirical estimations of rini, the offset of our 
Hl-masses by a factor ^ 2 can be readily accounted for 
by artificially decreasing the flux limit of the simulated 
survey by a factor 2. For typical Hl-surveys in the red- 

^ In the simulati on, both masses and lengths scale as 
HSpringel et al.l[2005h . 
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Fig. 7. — Comparison of the dN / dz-\>\ots for different emission lines observed with an identical peak flux density limit of 1 /iJy. The left 
panel shows our simulation presented in this paper. For comparison, the right panel represents the case of a simulation with no galaxy 
evolution, as obtained by using only the local galaxy simulation-box for the construction of the mock observing cone. The line types are 
as explained in Fig. (5] and colors have been used to distinguish HI (blue) from CO (red). All simulated functions (solid and dotted) use a 
redshift bin size of Az = 0.1 and a sky field of 4 X 4 deg^, corresponding to a maximal redshift Zmax = 10. 



shift range z — 0.5 — 10, this would increase the number 
of detectable sources by a factor 2 — 4. 

We have limited our predictions for HI to Hl-emission 
from galaxies. However, in the EoR, the IGM was not 
completely ionized and therefore ac ted as an additio nal 
source of Hl-emission or -absorption (jlliev et al.l2002[ ). It 
may therefore be necessary to analyze the implications 
of inte rgalactic HI on the detectability of galactic HI at 
z > 6 (jBecker et al.|[200l[l . On a theoretical level, such 
an analysis could result from combining the simulation 
presented in this paper with a simula tion of the EoR (e.g. 
iBaek et all [20091: iSantos et al.|[2008l ). 

4.3. Limitations specific to the CO-lines 

The main discussion of these limitations is given in 
lObreschkow fc Rawlingsl (|2009b( ). where we introduced 
our model for the conversion between H2 and CO. The 
most serious sources of uncertainty appear to be the heat- 
ing of molecular gas by SBs and AGNs, the overlap of 
molecular clouds at high redshift, and the possible pres- 
ence of nuclear molecular disks in high-redshift galax- 
ies. By contrast, the often discussed effects of the CMB 
and the cosmic evolution of the metallicity seem rela- 
tively well understood today. Overall, the uncertainty 
in the predicted CO-luminosities increases with redshift 
and with the J-level of the CO-transition. 

The highest uncertainties, i.e. those for the higher or- 
der CO-lines at high redshift, can be close to a factor 
10. The dA'^/dz-functions in this regime are therefore ex- 
pected to deviate significantly from our predictions. Such 
deviations will uncover much of the physics of CO-line 
emission. In fact, in Obreschkow & Rawlings (2009b) 
we have explained in detail how different deviations of 
the CO-luminosity functions from our predictions can be 
translated into physical interpretations. 

4.4. Is the simulation large enough to probe cosmic 
structure? 

The largest coherence-scale of our sky simulation is 
defined by the size of the periodic simulation box of the 
underlying dark matter simulation (Millennium Simula- 
tion, ISpringcl et al...2005<) . The side length of this box is 



Sbox = 500 Mpc, which sets the smallest extractable 
wave number to fc = 27r/sbox ~ 0.013 /i. This value is 
comparable to the wave nu mber of the first pea k in the 
CDM power spectrum (e.g. iSpringel et ani2005 ?). There- 
fore, the presented simulation allows us to study the 
power spectrum of HI- and CO-lines and to extract the 
baryon acoustic oscillations (BAOs); however the posi- 
tion and the amplitude of the first peak of the BAOs will 
be very poorly constrained. 

By contrast, the SKA will have the potential to im- 
prove on present measurements of the baryonic power 
spectrum by at least an order of magnitude in ampli- 
tude, and it will detect power in spatial frequencies far 
below the first acoustic peak. Such a detection could set a 
primordial constraint on cosmological para meters, espe- 
cially on the equation of state of dark energy (jBlake et al.l 
[2OOI lAbdalla et al.l 120091 ) . Therefore, a' simulation of 
such a detection is regarded as a necessary step in design- 
ing the SKA. Yet, this requirement represents a major 
challenge since no current simulation of cosmic structure 
is large enough to accurately follow the largest acous- 
tic oscillations, while simultaneously resolving structures 
small enough to allow the assembly of typical galaxies. 

A circumvention of this numerical predicament could 
result from merging two simulations w ith different 
length-scales (see e.g. lAngulo et al.l 120081 ). For exam- 
ple, we could adopt the Horizon-47r dark matt er sim- 
ulation (|Prunet et all 120081 : iTbvssier et al.l 120081 ) . which 
has a giant box side length of Sbox = 2/i~^Gpc, yet 
10-times less mass resolution than the Millennium Sim- 
ulation. Each dark matter halo of the Horizon-47r Sim- 
ulation could then be populated with the resolved dark 
matter substructure and the galaxies contained in com- 
parable haloes of the Millennium Simulation. 

5. CONCLUSION 

With this paper we release a simulation of the Hl- 
emission line and the first ten ^^C^^O-emission lines 
of galaxies in a sky field with a comoving diameter of 
500 h~^ Mpc. The actual field-of-view depends on the 
(user-defined) maximal redshift Zmax according to the 
relation displayed in Fig. [2] (see also Eq. [5]) . This Simula- 
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tion represents the first quantitative attempt to compare 
the detectability of HI and CO at high redshift. De- 
spite the limitations and uncertainties of this simulation 
(Section |4]) , its underlying galaxy simulation is nonethe- 
less consisten t with currently avai lable observations 
Csee lObreschkow et alJ l2 009at lObreschkow fc Rawlingsl 
I2009bt lObreschkow et a l. 2009b). 

While this paper focussed on the simulation techniques 
and directly accessible results, the list of possible appli- 
cations of the presented simulation is extensive. Some 
examples are: 

1. a dA^/d2;-analysis for particular surveys with the 
SKA, the LMT, and ALMA or their pathfinders; 

2. a combined study of Hl-emission from galaxies 
and Hl-emission fr om the IGM during the EoR 
(jSantos et al.ll2008l) to quantify confusion issues; 

3. an optimization of the survey time allocated to dif- 
ferent ALMA bands based on the CO-line ratios 
predicted at various redshifts; 

4. first quantitative predictions of the SKA's and 
ALMA's abilities to probe the galaxy power spec- 
trum; 

5. predictions of the absorption signatures of HI and 
CO against distant QSOs; 



6. a study of line stacking experiments at redshifts 
where the detection of individual galaxies becomes 
impossible. 



Such predictions can assist the design and optimized 
use of telescopes, such as the LMT, the SKA, ALMA, 
and their pathfinders. Moreover, in light of forthcoming 
observations with those future instruments, the predic- 
tions made prior to these observations are the safest and 
perhaps the only way to test the predictive power of our 
current theories. This feature accentuates the necessity 
of extensive predictions, whether they will be verified or 
not by the empirical data. 



This effort /activity is supported by the European 
Community Framework Programme 6, Square Kilome- 
tre Array Design Studies (SKADS), contract no 011938. 
The Millennium Simulation databases and the web appli- 
cation providing online access to them were constructed 
as part of the activities of the German Astrophysical Vir- 
tual Observatory. 
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APPENDIX 

A. HOW TO RECOVER LINE PROFILES FROM CATALOG PROPERTIES 

In Section 12. 4[ the normalized profiles of the HI- and CO-emission lines have been parameterized using the five 
parameters w°ll^, and W2o^ (see Fig. [8]). From these parameters, the original normalized velocity 

profiles "^(V) can be approximately recovered using the analytic function 



-0.5 



if 11^1 > <ik/2, 



^obs 



(Al) 



if 11^1 < <lk/2 and > 
if 1^1 < <lk/2 and = 

where fc ,;, i € {1, 5}, are free parameters. Eq. (|Aip combines the the functional form of Eq. (42) in lObreschkow et al.l 
(|2009af ) for the center of the emission line with exponential tails. The five parameters ki can be inferred from the 
parameters ^q^^, ^I'max; ""^peik' '"^50'': and 'W2o^. The analytic solution is 



ki- 

k2- 



-0.693-2.322 



0.842 



ln(it;: 

„,,obs 



obs 
20 



W. 



obs \ 
peak/ 



ln(w; 



iObs 
50 



w 



,obs \ 
peak/ 



peak 



4 1^0 



obs 
ma: 

3hs 



if > vE-g^ 



if "^"^ 

max 



(A2) 
(A3) 

(A4) 
(A5) 
(A6) 



Fig. [8] compares a simulated normalized Hl-emission line with the emission line recovered using Eqs. (jAlHAGI 




Fig. 8. — Comparison of a simulated normalized emission line '^(V) (solid line) with the emission line 'l'approx(V) (dashed line), recovered 



from the five parameters , 'I'n 



peak ' 



^, and w'. 



obs 
20 



We emphasize that the approximation of Eq. (jAll) does not conserve the normalization condition '^{V) — 1- 
However, the integral ^'approx(^) normally differs from unity by no more than 10%. Hence, a good approximation 
of the flux proflle can be obtained by multiplying 4'approx(^) with the velocity-integrated line fiux of the respective 
emission line as given in the simulation catalog. 

We note that some emissio n line s, especially those of galaxies seen face-on, peak at the line center. These lines have 
^pcak ~ 0, and therefore Eq. (|Aip reduces to the exponential tails. 
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B. ONLINE-ACCESS TO THE SKY SIMULATION 

One particular realization of the mock observing cone, i.e. one choice of random symmetry operations for the 
replicated simulation boxes (see Section can be access on-line via http://s-cubed.physics.ox.ac.uk/ (go to "S'^- 
SAX-Sky"). Each galaxy in the virtual observing cone is specified by a list of properties, including its position, its 
attributes for the HI- and CO-emission lines, as well as its intrinsic properties of the DeLucia-catalog, such as optical 
magnitudes, masses, star formation rates, clustering properties, or merger histories. 

The accessible database contains two subsets, a full sky simulation (^ 2.8 • 10* galaxies), associated with the 
Millennium Simulation (sbox = 500/i~^Mpc), and a small sky simulation (~ 4.5 • 10^ galaxies), associated with the 
Milli-Millennium Simulation (sbox = 62.5 Mpc). The solid angle subtended by the small sky simulation is 64-times 
smaller than of the full sky simulation. However, the small simulation can be useful for testing purposes, since it can 
be accessed and post-processed 10-100 times faster. 

The maximal opening angle (p of the sky field depends on the maximal comoving distance I3c,max (or the maximal 
redshift Zmax) via Eq. [5] (see also Fig. [2]) . The user must be aware that there are no galaxies outside this maximal 
opening angle, i.e. galaxies only exist, where |RA| and |Dec| are smaller than (p/2. 

The database can be queried using the structured query language (SQL) interface. The latter not only allows 
the user to download a particular galaxy sample, but it also offers ways to directly calculate dA/dz-functions and 
luminosity-functions or to retrieve all the galaxies of a particular cluster. Samples of such advanced queries are given 
on the web-page. 
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C. ILLUSTRATION OF A SKY FIELD IN DIFFERENT REDSHIFT BINS 

Fig. |9] shows the HI and CO of the galaxies in a mock sky field of 3 x 1 arcmin^ at three different redshifts. Each 
redshift slice has the same comoving thickness of 240 Mpc, such that the number of galaxies is proportional to the 
comoving space density of galaxies. Note, however, that the flux scales differ between the three panels of Fig. [9l The 
main features of the galaxies in Fig. [HI are discussed in Section 13.11 




Fig. 9. — Simulated sky field covering 3x1 arcmin^. The three panels correspond to three different redshift slices, with an identical 
comoving depth of 240 Mpc. The coloring is identical to Fig. [l] but the flux scales are 10-times smaller at ^ = 3 and 100-times smaller at 
z = 6. 
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D. ANALYTIC FITS FOR dA^/dz-FUNCTIONS 

Table [T] lists the values of the parameters ci, C2, and C3 for the analytic dTV/dz-fit of Eq. (fT5|) . The parameters are 
given for both peak flux density limited and integrated flux limited surveys, although the functions have only been 
displayed for the peak flux density limited case. Since apparent line widths are typically of order lOOkms"^, surveys 
limited by a flux density sum and those limited by a velocity-integrated flux = siim x lOOkms"^ have similar 
dN/ dz-functions. 

The sensitivity of a survey may be high enough that galaxies are detected in regions of the galaxy MF, where the 
simulation is incomplete. In fact, for each emission line and each sensitivity limit, there is a critical redshift Zc, below 
which sources of the incomplete part of the simulated galaxy MF will be detected at a sufficiently high rate that 
the number of real detections will signiflcantly exceed the number of simulated detections. Therefore, the simulated 
diV/dz-functions at z < Zc must be considered as lower limits. As in Section [321 we shall deflne Zc as the redshift, 
below which the fraction of sources with hydrogen masses below 10^ M0 (completeness limit) is larger than 1%, noting 
that the fraction of actually missing sources may be much larger due to incompleteness (see Section [3.2p . 



Parameters 



Limiting 

lQ-8 



peak 

10"^ 



flux density 



[Jy] 



10" 



Limiting 

10"^ 10"^ 



integrated flux 

10"* 10"^ 



[Jy km ! 
10"2 



10" 



HI 



6.55 
2.54 
1.42 
2.3 



6.87 
2.85 
2.17 
1.1 



6.73 
2.32 
3.09 
0.5 



5.75 
1.14 
3.95 
0.2 



4.56 
0.43 
6.86 
0.1 



6.62 
2.64 
35.49 
0.1 



6.61 
2.66 
1.51 
2.1 



6.94 
2.91 
2.34 
0.9 



6.53 
2.09 
2.86 
0.4 



5.79 
1.33 
3.92 
0.1 



5.21 
1.23 
8.38 
0.1 



6.33 
2.60 
30.52 
0.1 



CO(l-O) 



6.32 
2.16 
1.11 
6.0 



6.40 
2.29 
1.27 
5.1 



6.37 
2.05 
1.46 
2.2 



5.92 
1.33 
1.55 
0.9 



5.06 
0.54 
1.72 
0.3 



4.15 
0.31 
3.46 
0.1 



6.31 
2.14 
1.11 
6.0 



6.37 
2.23 
1.23 
5.1 



6.30 
1.96 
1.38 
2.2 



5.91 
1.43 
1.47 
0.9 



5.27 
0.95 
1.69 
0.4 



5.17 
1.37 
5.19 
0.1 



CO(2-l) 



6.29 
2.10 
1.07 
5.8 



6.34 
2.20 
1.15 
6.0 



6.42 
2.27 
1.34 
3.8 



6.25 
1.81 
1.51 
1.5 



5.64 
1.07 
1.58 
0.6 



4.64 
0.35 
1.94 
0.2 



6.29 
2.09 
1.06 
5.8 



6.33 
2.18 
1.14 
6.0 



6.38 
2.18 
1.28 
3.8 



6.21 
1.84 
1.46 
1.5 



5.67 
1.23 
1.49 
0.6 



4.87 
0.70 
1.82 
0.2 



CO(3-2) 



6.28 
2.08 
1.06 
5.8 



6.32 
2.16 
1.12 
6.0 



6.40 
2.27 
1.27 
4.9 



6.32 
1.95 
1.44 
2.0 



5.79 
1.22 
1.48 
0.7 



4.86 
0.47 
1.58 
0.3 



6.28 
2.08 
1.06 
5.8 



6.31 
2.14 
1.11 
6.0 



6.37 
2.21 
1.23 
4.7 



6.25 
1.89 
1.36 
2.0 



5.81 
1.36 
1.43 
0.8 



5.09 
0.87 
1.56 
0.3 



CO (4-3) 



6.28 
2.08 
1.06 
5.8 



6.31 
2.13 
1.10 
6.0 



6.37 
2.22 
1.21 
5.2 



6.26 
1.87 
1.32 
2.1 



5.74 
1.18 
1.36 
0.7 



4.84 
0.52 
1.38 
0.3 



6.28 
2.08 
1.06 
5.8 



6.30 
2.12 
1.09 
6.0 



6.34 
2.17 
1.18 
5.2 



6.19 
1.79 
1.25 
2.0 



5.77 
1.34 
1.33 
0.8 



5.03 
0.86 
1.36 
0.3 



CO (5-4) 



6.28 
2.07 
1.05 
5.8 



6.30 
2.11 
1.09 
6.0 



6.32 
2.13 
1.17 
5.1 



6.05 
1.56 
1.15 
1.7 



5.45 
0.96 
1.14 
0.6 



4.53 
0.54 
1.14 
0.2 



6.28 
2.07 
1.05 
5.8 



6.29 
2.11 
1.08 
6.0 



6.30 
2.09 
1.15 
4.7 



6.00 
1.52 
1.09 
1.6 



5.53 
1.20 
1.17 
0.6 



4.74 
0.83 
1.15 
0.2 



CO (6-5) 



6.28 
2.07 
1.05 
5.8 



6.29 
2.09 
1.09 
6.0 



6.17 
1.81 
1.07 
6.2 



5.62 
1.04 
0.89 
0.9 



4.95 
0.78 
0.92 
0.3 



4.17 

0.93 
1.08 
0.1 



6.28 
2.07 
1.05 
5.8 



6.28 
2.07 
1.08 
6.0 



6.13 
1.76 
1.03 
6.0 



5.62 
1.11 
0.87 
0.9 



5.08 
1.02 
0.97 
0.4 



4.35 
1.13 
1.06 
0.1 



CO (7-6) 



6.27 
2.06 
1.05 
6.0 



6.22 
1.96 
1.06 
5.2 



5.72 
1.14 
0.80 
1.3 



5.01 
0.69 
0.65 
0.4 



4.45 
1.13 
0.88 
0.2 



3.91 
1.58 
1.24 
0.1 



6.27 
2.05 
1.05 
6.0 



6.20 
1.93 
1.04 
5.2 



5.70 
1.18 
0.79 
1.2 



5.08 
0.80 
0.66 
0.4 



4.55 
1.18 
0.86 
0.1 



4.06 
1.72 
1.14 
0.1 



CO (8-7) 



6.22 
1.97 
1.04 
6.0 



5.76 
1.25 
0.78 
6.6 



5.00 
0.66 
0.51 
0.5 



4.42 
1.04 
0.62 
0.2 



4.13 
1.85 
1.03 
0.1 



3.68 
2.22 
1.45 
0.1 



6.20 
1.95 
1.02 
6.0 



5.74 
1.27 
0.78 
1.3 



5.08 
0.80 
0.55 
0.5 



4.48 
0.99 
0.58 
0.2 



4.19 
1.83 
0.97 
0.1 



3.78 
2.06 
1.19 
0.1 



CO (9-8) 



5.70 
1.24 
0.73 
6.2 



4.95 
0.74 
0.46 
0.4 



4.35 
1.05 
0.48 
0.1 



4.10 
1.81 
0.81 
0.1 



3.84 
2.31 
1.14 
0.1 



3.33 
2.17 
1.43 
0.1 



5.69 
1.27 
0.73 
1.2 



5.03 
0.87 
0.51 
0.4 



4.43 
1.04 
0.48 
0.1 



4.12 
1.70 
0.73 
0.1 



3.93 
2.20 
1.04 
0.1 



3.43 
2.27 
1.21 
0.1 



CO (10-9) 



4.81 
0.86 
0.42 
0.3 



4.29 
1.29 
0.48 
0.1 



4.03 
1.80 
0.66 
0.1 



3.88 
2.29 
0.97 
0.1 



3.52 
2.27 
1.11 
0.1 



2.95 
2.23 
1.47 
0.1 



4.89 
0.95 
0.47 
0.3 



4.35 
1.26 
0.48 
0.1 



4.07 
1.76 
0.64 
0.1 



3.91 
2.21 
0.89 
0.1 



3.65 
2.35 
1.06 
0.1 



3.00 
2.03 
1.06 
0.1 



TABLE 1 

Parameters for the analytic fit formula of Eq. (|15|l for dN/dz peak flu.x density limited and integrated flux limited 

SURVEYS. 



